% Calculating CIELAB values for reflectances and cmf weighted with D50.
clear all
reflectance = xlsread('Data_color__samples.xls','B4:AF27');

cmf = xlsread('Data_color__samples.xls','Weights','B5:D35'); % Actually source*cmf (d50 at B5:D35)

% Normalize reflectance to normal-people numbers!
reflectance = reflectance ./ 100;

% Calculate source*reflectance*cmf for each sample
xyz_d50 = reflectance * cmf;

% Calculate tristimulus for the white point automagically!
wp = sum(cmf, 1);

for i = 1:size(xyz_d50,1)
    % calculating lab for every xyz value.
    % observer is always accomodated to the source
 cielab_d50(i,:) = xyz2lab(xyz_d50(i,:),wp);
end

cmf = xlsread('Data_color__samples.xls','Weights','G5:I35'); % Actually source*cmf (d50 at B5:D35)

% Calculate source*reflectance*cmf for each sample
xyz_d65 = reflectance * cmf;

% Calculate tristimulus for the white point automagically!
wp = sum(cmf, 1);

for i = 1:size(xyz_d65,1)
    % calculating lab for every xyz value.
    % observer is always accomodated to the source
 cielab_d65(i,:) = xyz2lab(xyz_d65(i,:),wp);
end

cmf = xlsread('Data_color__samples.xls','Weights','L5:N35'); % Actually source*cmf (d50 at B5:D35)

% Calculate source*reflectance*cmf for each sample
xyz_a = reflectance * cmf;

% Calculate tristimulus for the white point automagically!
wp = sum(cmf, 1);

for i = 1:size(xyz_a,1)
    % calculating lab for every xyz value.
    % observer is always accomodated to the source
 cielab_a(i,:) = xyz2lab(xyz_a(i,:),wp);
end

cmf = xlsread('Data_color__samples.xls','Weights','Q5:S35'); % Actually source*cmf (d50 at B5:D35)

% Calculate source*reflectance*cmf for each sample
xyz_d75 = reflectance * cmf;

% Calculate tristimulus for the white point automagically!
wp = sum(cmf, 1);

for i = 1:size(xyz_d75,1)
    % calculating lab for every xyz value.
    % observer is always accomodated to the source
 cielab_d75(i,:) = xyz2lab(xyz_d75(i,:),wp);
end

sRGBMat = [3.2404542 -1.5371385 -0.4985314;...
        -0.9692660  1.8760108  0.0415560;...
        0.0556434 -0.2040259  1.0572252];
sRGB_d50 = sRGBMat * xyz_d50'/100.*255;
sRGB_d50(sRGB_d50 < 0) = 0;
sRGB_d65 = sRGBMat * xyz_d65'/100*255;
sRGB_d65(sRGB_d65 < 0) = 0;
sRGB_d75 = sRGBMat * xyz_d75'/100*255;
sRGB_d75(sRGB_d75 < 0) = 0;
sRGB_a = sRGBMat * xyz_a'/100*255;
sRGB_a(sRGB_a < 0) = 0;



%%
close all;

figure;
hold on;
grid on;
box on;
plot(cielab_d50(:,2), cielab_d50(:,3), '*b');
plot(cielab_a(:,2), cielab_a(:,3), '*r');
plot([-100 100],[0 0],'k','LineWidth',2);
plot([0 0],[-100 100],'k','LineWidth',2);
arrow(cielab_d50(:, 2:3), cielab_a(:, 2:3), 'Length', 8, 'TipAngle', 30);
axis equal;
axis([-60 60 -80 100]);
set(gca, 'XTick', -60:20:60);
set(gca, 'YTick', -80:20:100);
xlabel('a*');
ylabel('b*');
title('Chroma Differences Between Illuminants D50 and A');
legend('D50', 'A', 'Location', 'NorthWest');

print('-depsc','d50_vs_a')

figure;
hold on;
grid on;
box on;
plot(cielab_d75(:,2), cielab_d75(:,3), '*b');
plot(cielab_a(:,2), cielab_a(:,3), '*r');
plot([-100 100],[0 0],'k','LineWidth',2);
plot([0 0],[-100 100],'k','LineWidth',2);
arrow(cielab_d75(:, 2:3), cielab_a(:, 2:3), 'Length', 8, 'TipAngle', 30);
axis equal;
axis([-60 60 -80 100]);
set(gca, 'XTick', -60:20:60);
set(gca, 'YTick', -80:20:100);
xlabel('a*');
ylabel('b*');
title('Chroma Differences Between Illuminants D75 and A');
legend('D75', 'A', 'Location', 'NorthWest');

print('-depsc','d75_vs_a')

figure;
hold on;
grid on;
box on;
plot(cielab_d75(:,2), cielab_d75(:,3), '*b');
plot(cielab_d50(:,2), cielab_d50(:,3), '*r');
plot([-100 100],[0 0],'k','LineWidth',2);
plot([0 0],[-100 100],'k','LineWidth',2);
arrow(cielab_d75(:, 2:3), cielab_d50(:, 2:3), 'Length', 8, 'TipAngle', 30);
axis equal;
axis([-60 60 -80 100]);
set(gca, 'XTick', -60:20:60);
set(gca, 'YTick', -80:20:100);
xlabel('a*');
ylabel('b*');
title('Chroma Differences Between Illuminants D75 and D50');
legend('D75', 'D50', 'Location', 'NorthWest');

print('-depsc','d75_vs_d50')

patchnames = {'Dark Skin', 'Light Skin', 'Blue Sky', 'Foliage', 'Blue Flower', ...
    'Bluish Green', 'Orange', 'Purplish Blue', 'Moderate Red', 'Purple', ...
    'Yellow Green', 'Orange Yellow', 'Blue', 'Green', 'Red', 'Yellow', 'Magenta',...
    'Cyan', 'Neutral 9.5', 'Neutral 8', 'Neutral 6.5', 'Neutral 5', 'Neutral 3.5', ...
    'Neutral 2'};


% Just trust me. This works.
file_h = fopen('data_d50_a.txt', 'w');
fprintf(file_h, '\\begin{tabular}{|c|c|c|c|c|c|c|c|c|}\n');
fprintf(file_h, '\\hline\n');
fprintf(file_h, ' & \\multicolumn{3}{|c|}{$D_{50}$} & \\multicolumn{3}{|c|}{$A$} & \\multicolumn{2}{|c|}{Color Differences} \\\\ \\hline\n');
fprintf(file_h, 'Color Patches & L* & a* & b* & L* & a* & b* & $\\Delta E*_{ab}$ & $\\Delta E_{94}$ \\\\ \\hline\n');
for i = 1:24
    fprintf(file_h, '%s & %.2f & %.2f & %.2f & %.2f & %.2f & %.2f & %.2f & %.2f \\\\ \\hline\n',...
        patchnames{i}, cielab_d50(i,1), cielab_d50(i,2), cielab_d50(i,3),...
        cielab_a(i,1), cielab_a(i,2), cielab_a(i,3),...
        getDelab(cielab_d50(i,:), cielab_a(i,:)),...
        getDel94(cielab_d50(i,:), cielab_a(i,:)));
end
fprintf(file_h, '\\end{tabular}');
fclose(file_h);

% Just trust me. This works.
file_h = fopen('data_d75_a.txt', 'w');
fprintf(file_h, '\\begin{tabular}{|c|c|c|c|c|c|c|c|c|}\n');
fprintf(file_h, '\\hline\n');
fprintf(file_h, ' & \\multicolumn{3}{|c|}{$D_{75}$} & \\multicolumn{3}{|c|}{$A$} & \\multicolumn{2}{|c|}{Color Differences} \\\\ \\hline\n');
fprintf(file_h, 'Color Patches & L* & a* & b* & L* & a* & b* & $\\Delta E*_{ab}$ & $\\Delta E_{94}$ \\\\ \\hline\n');
for i = 1:24
    fprintf(file_h, '%s & %.2f & %.2f & %.2f & %.2f & %.2f & %.2f & %.2f & %.2f \\\\ \\hline\n',...
        patchnames{i}, cielab_d75(i,1), cielab_d75(i,2), cielab_d75(i,3),...
        cielab_a(i,1), cielab_a(i,2), cielab_a(i,3),...
        getDelab(cielab_d75(i,:), cielab_a(i,:)),...
        getDel94(cielab_d75(i,:), cielab_a(i,:)));
end
fprintf(file_h, '\\end{tabular}');
fclose(file_h);

% Just trust me. This works.
file_h = fopen('data_d75_d50.txt', 'w');
fprintf(file_h, '\\begin{tabular}{|c|c|c|c|c|c|c|c|c|}\n');
fprintf(file_h, '\\hline\n');
fprintf(file_h, ' & \\multicolumn{3}{|c|}{$D_{75}$} & \\multicolumn{3}{|c|}{$D_{50}$} & \\multicolumn{2}{|c|}{Color Differences} \\\\ \\hline\n');
fprintf(file_h, 'Color Patches & L* & a* & b* & L* & a* & b* & $\\Delta E*_{ab}$ & $\\Delta E_{94}$ \\\\ \\hline\n');
for i = 1:24
    fprintf(file_h, '%s & %.2f & %.2f & %.2f & %.2f & %.2f & %.2f & %.2f & %.2f \\\\ \\hline\n',...
        patchnames{i}, cielab_d75(i,1), cielab_d75(i,2), cielab_d75(i,3),...
        cielab_a(i,1), cielab_d50(i,2), cielab_d50(i,3),...
        getDelab(cielab_d75(i,:), cielab_d50(i,:)),...
        getDel94(cielab_d75(i,:), cielab_d50(i,:)));
end
fprintf(file_h, '\\end{tabular}');
fclose(file_h);


